
float maxDist = 10;
float myCir1P[] = new float[]{0,0,0};
float myCir2P[];
float myR1=50,myR2=50;

void drawIntersect(float p1[],float p2[],float r1,float r2) {
  int i=0,j=0,k=0;
  float[] midPoint = new float[]{(p1[0]+p2[0])/2,(p1[1]+p2[1])/2,(p1[2]+p2[2])/2};
  float[] r1TOr2V = new float[]{p2[0]-p1[0],p2[1]-p1[1],p2[2]-p1[2]};
  float distance = DogaMakeUnitVector(r1TOr2V);
  float rN=0;
//  distance = sqrt(distance);
  if (distance>(r1+r2)/2) return;
  rN = -distance+r1+r2;
  rN = rN/(r1+r2);
  rN = rN*sqrt(r1*r1+r2*r2);
  //draw a circle with radius rN at the point midPoint
//  fill(255,0,255);
  ellipse(midPoint[0],midPoint[1],rN,rN);
}
void setup(){
  size(400,400);
  myCir2P = new float[3];
  myCir2P[0] = width/2;
  myCir2P[1] = height/2;
  myCir2P[2] = 0;
  ellipseMode(CENTER);
  noStroke();
}

void draw() {
  background(0,255,0);
  fill(0,0,255);
  ellipse(myCir1P[0],myCir1P[1],myR1,myR1);
//  fill(255,0,0);
  ellipse(myCir2P[0],myCir2P[1],myR2,myR2);
  drawIntersect(myCir1P,myCir2P,myR1,myR2);
}

void mousePressed() {
  myCir1P[0] = mouseX;
  myCir1P[1] = mouseY; 
}
float sinus (float degree) {return (float)sin(6.28*(float)degree/360);}
float cosinus (float degree) {return (float)cos(6.28*(float)degree/360);}
float rad2deg (float radian) {
  while (radian<0) radian+=6.28;
  while (radian>6.28) radian-=6.28;
  return 360*(radian/6.28);
}
float absD(float x) {
  if (x>=0) return x;
  return -x;
}
/***********************************************************/
/** B A S I C  V E C T O R  &  M A T R I X  F U N C T I O N S **/
/***************************************************************/
void DogaPrintMatrixCol(String matrixName,float[][] Matrix,int row,int col/*PRINTS IN COL ORDER*/) {
  println("Printing the matrix = "+matrixName+" in coloumn based order");
  for (int i=0;i<row;i++) {
	for (int j=0;j<col;j++) {
		print("["+Matrix[i+j*row]+"],");
	}
	println();
  }
  println();
  return;
}
void DogaPrintMatrixRow(String matrixName,float[][] Matrix,int row,int col/*PRINTS IN COL ORDER*/) {
  println("Printing the matrix = "+matrixName+" in row based order");
  for (int i=0;i<row;i++) {
	for (int j=0;j<col;j++) {
		print("["+Matrix[j+i*col]+"],");
	}
	println();
  }
  println();
  return;
}
void DogaCopyMatrix(float[] from,float[] to,int n){
  int i=0;
  for (i=0;i<n;i++) to[i]=from[i];
  return;
}
void DogaCopyMatrix(float[][] from,float[][] to,int row,int col){
  int i=0,j=0;
  for (i=0;i<row;i++)
    for(j=0;j<col;j++)
      to[i][j]=from[i][j];
  return;
}
float DogaDotProduct(float V1[],float V2[]) {return V1[0]*V2[0]+V1[1]*V2[1]+V1[2]*V2[2];}
float[] DogaCrossProduct(float V1[],float V2[]) {
  float m;
  float[] Vn=new float[4];
  Vn[0] = V1[1]*V2[2] - V2[1]*V1[2];
  Vn[1] = V2[0]*V1[2] - V1[0]*V2[2];
  Vn[2] = V1[0]*V2[1] - V2[0]*V1[1];
  Vn[3] = 1;
  m = sqrt(DogaDotProduct(Vn,Vn));
  if (m==0) {
    println("**BE CAREFUL THE LINES ARE PARALLEL TO EACH OTHER***");
    for (int i=0;i<3;i++) Vn[i] = 0;
  }
  else for (int i=0;i<3;i++) Vn[i] = Vn[i]/m;
  return Vn;
}
float DogaCrossProduct(float V1[],float V2[],float ApplySolution[]) {
  float m;
  float[] Vn=new float[4];
  Vn[0] = V1[1]*V2[2] - V2[1]*V1[2];
  Vn[1] = V2[0]*V1[2] - V1[0]*V2[2];
  Vn[2] = V1[0]*V2[1] - V2[0]*V1[1];
  Vn[3] = 1;
  m = sqrt(DogaDotProduct(Vn,Vn));
  if (m==0) {
    println("**BE CAREFUL THE LINES ARE PARALLEL TO EACH OTHER***");
    for (int i=0;i<3;i++) ApplySolution[i] = 0;
  }
  else for (int i=0;i<3;i++) ApplySolution[i] = Vn[i]/m;
  ApplySolution[3]=1;
  return m;
}
float DogaMakeUnitVector(float V[]) {
  float w=0;
  w = sqrt(DogaDotProduct(V,V));
  if (w!=0) {
    V[0] = V[0]/w;
    V[1] = V[1]/w;
    V[2] = V[2]/w;
  }
  else {
    V[0] = 0;
    V[1] = 0;
    V[2] = 0;
  }
  return w;
}
float[] DogaReturnUnitVector(float V[]) {
  float w=0,retV[]=new float[4];
  w = sqrt(DogaDotProduct(V,V));
  if (w!=0) {
    retV[0] = V[0]/w;
    retV[1] = V[1]/w;
    retV[2] = V[2]/w;
  }
  else {
    retV[0] = 0;
    retV[1] = 0;
    retV[2] = 0;
  }
  retV[3] = 1;
  return retV;
}
void DogaLoadIdentity(float[][] Matrix,int row,int col) {
  int i=0,j=0;
  for (i=0;i<row;i++)
    for (j=0;j<col;j++)
      if (i==j) Matrix[i][j]=1; else Matrix [i][j]=0;
}

float[][] DogaMatrixMultCol(float[][] matA,float[][] matB,int A,int B,int C) {
  float sol[][]=new float[A][C];
  int i,j,k;
  for (i=0;i<A;i++) {
    for (j=0;j<C;j++) {
      sol[i][j] = 0;
      for (k=0;k<B;k++)
  	sol[i][j] += matA[i][k]*matB[k][j];
    }
  }
  return sol;
}

void DogaMatrixMultCol(float[][] matA,float[][] matB,int A,int B,int C,float[][] ApplySolution) {
  float sol[][]=new float[A][C];
  int i,j,k;
  for (i=0;i<A;i++) {
    for (j=0;j<C;j++) {
      sol[i][j] = 0;
      for (k=0;k<B;k++)
  	sol[i][j] += matA[i][k]*matB[k][j];
    }
  }
  for (i=0;i<A;i++)
    for (j=0;j<C;j++)
       ApplySolution[i][j]=sol[i][j]; 
  return;
}
float[][] DogaMatrixMultRow(float[][] matA,float matB[][],int A,int B,int C) {
  float sol[][]=new float[A][C];
  int i,j,k;
  // matA will be AxB matrix ; matB will be BxC matrix ; solution will be AxC matrix
  for (i=0;i<A;i++)
    for (j=0;j<C;j++)
      sol[i][j] = 0;
  for (i=0;i<A;i++) {//A rows(points) of matA
    for (j=0;j<B;j++) {//B coloumns of matA and rows of matB
      for (k=0;k<C;k++) //C cols of matB
	sol[i][j] += matA[i][k]*matB[j][k];
    }
  }
  return sol;
}
void DogaMatrixMultRow(float[][] matA,float matB[][],int A,int B,int C,float[][] ApplySolution) {
  float sol[][]=new float[A][C];
  int i,j,k;
  // matA will be AxB matrix ; matB will be BxC matrix ; solution will be AxC matrix
  for (i=0;i<A;i++)
    for (j=0;j<C;j++)
      sol[i][j] = 0;
  for (i=0;i<A;i++) {//A rows(points) of matA
    for (j=0;j<B;j++) {//B coloumns of matA and rows of matB
      for (k=0;k<C;k++) //C cols of matB
	sol[i][j] += matA[i][k]*matB[j][k];
    }
  }
  for (i=0;i<A;i++)
   for (j=0;j<C;j++)
     ApplySolution[i][j] = sol[i][j];
  return;
}

/*
  float r1TOr2DIFX=p2[0]-p1[0];
  float r1TOr2DIFY=p2[1]-p1[1];
  float r1TOr2DIFZ=p2[2]-p1[2];
  float r1TOr2V[] = float[]{p2[0]-p1[0],p2[1]-p1[1],p2[2]-p1[2]};
  float r2TOr1DIFX=-r1TOr2DIFX;
  float r2TOr1DIFY=-r1TOr2DIFY;
  float r2TOr1DIFZ=-r1TOr2DIFZ;
  float r2TOr1V[] = float[]{p1[0]-p2[0],p1[1]-p2[1],p1[2]-p2[2]};
  float rMin = 0;
  DogaMakeUnitVector(r1TOr2V);
  DogaMakeUnitVector(r2TOr1V);  
*/
